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Development  of  a  Versatile  Physics-Based 
Finite-Element  Model  of  an  AlGaN/GaN  HEMT 
Capable  of  Accommodating  Process  and  Epitaxy 
Variations  and  Calibrated  Using 
Multiple  DC  Parameters 

Eric  R.  Heller,  Rama  Vetury,  and  Daniel  S.  Green,  Member,  IEEE 


Abstract — We  present  a  physics-based  finite-element  model  of 
operation  of  an  AlGaN/GaN  HEMT  with  device  geometry  in¬ 
puts  taken  from  transmission  electron  microscope  cross  sections 
and  calibrated  by  comparison  with  measured  electrical  data 
comprising  standard  field-effect  transistor  metrics  and  less  well- 
known  model  parameters.  A  variety  of  electrical  outputs  from  the 
model  are  compared  to  experiment,  and  the  level  of  agreement  is 
reported. 

Index  Terms — Device  model,  field-effect  transistor  (FET),  GaN, 
GaN/AlGaN,  high-electron  mobility  transistor  (HEMT),  model 
calibration,  model  characterization,  modulation-doped  field  effect 
transistor  (MODFET). 


I.  Introduction 

BECAUSE  of  the  high  breakdown  field,  high  band  gap, 
and  reasonably  high  mobility  of  GaN  and  the  AlGaN/ 
GaN  channel,  AlGaN/GaN  high-electron  mobility  transistors 
(HEMTs)  are  rapidly  being  commercialized  for  high-power  dc 
and  microwave  radio- frequency  devices.  In  fact,  because  of 
the  innate  promise  of  the  materials  system  and  the  progress 
already  made  toward  realizing  its  potential,  devices  are  being 
commercially  produced,  despite  the  presence  of  some  perfor¬ 
mance  issues  and  long-term  reliability  concerns  with  this  ma¬ 
terials  system  that  require  more  detailed  study.  Physics-based 
electrothermal  models  of  these  devices  have  been  constructed 
in  the  past  [l]-[4],  and  most  include  some  comparison  to  exper¬ 
imental  data  for  model  calibration  and/or  verification.  However, 
a  complete  and  convincing  device  model  should  be  able  to 
simultaneously  match  experimental  data  for  numerous  metrics 
and  should  realistically  simulate  the  thermal  characteristics  of 
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the  full  device  using  a  3-D  simulation  domain  that  includes  the 
whole  thermal  path. 

In  order  to  incorporate  all  of  the  physics  relevant  to  this  prob¬ 
lem,  we  have  constructed  a  finite-element  model  of  the  small- 
scale  physics  in  the  ISE  Sentaurus  Device  [5]  of  the  HEMT  and 
linked  this  to  a  large-scale  3-D  finite-element  thermal  model 
in  ANSYS  [6],  which  contains  the  larger  scale  features.  This 
allows  the  full  3-D  geometry  of  the  substrate  and  package  to  be 
included  in  the  model. 

II.  Model  and  Experimental  Details 

We  modeled  the  device  using  real  device  cross  sections,  with 
the  overall  shape  of  the  gate  extracted  from  TEM  cross  sections 
of  representative  devices  and  other  dimensions  extracted  from 
design  rules.  This  is  not  to  say  that  we  are  stating  that  this 
was  required  for  a  functioning  model;  we  just  mean  that  we 
used  the  most  accurate  data  available.  The  layer  stack  consisted 
of  an  AlGaN  barrier  structure  on  2.05-/im  GaN.  This  is  on  a 
100-nm  interfacial  layer  on  a  thinned  100-/xm  SiC  die.  The  part 
is  a  six-finger  370-/im-gate- width  AlGaN/GaN  HEMT  on  a  die 
measuring  800  /am  by  800  /am  mounted  on  12.7-/im  (1/2  mil) 
thermal  die  attach  (50.6  W/m/K)  in  a  CuW  package,  which  is 
held  at  a  fixed  temperature  after  1.54-mm  (60  mils)  thickness. 

Inputs  to  the  electrothermal  (small  scale)  portion  of  the 
device  model  include  temperature-dependent  thermal  conduc¬ 
tivities  [1],  [2],  [7]  and  heat  capacities;  temperature-  and 
mole-fraction-dependent  semiconductor  band  gaps,  and  ohmic 
contact  resistances;  dielectric  constants;  electron  effective 
masses  (for  gate  tunneling  leakage  and  density  of  states); 
hot  electron  relaxation  time  (Hydrodynamic  model)  [8];  and 
field-  and  temperature-dependent  electron  mobility  fitted  to 
Albrecht  et  al.  [9],  except  that  the  low  field  mobility  has  been 
set  to  1630  cm2/V/s  based  on  Hall  measurements  on  Van  der 
Pauw  test  patterns.  Additional  details  on  the  implementation  of 
these  inputs  can  be  found  in  our  earlier  publications  [1],  [2]  and 
the  user  guide  for  the  ISE  Sentaurus  Device  [5]. 

Gate  reverse  current  modeling  includes  the  effects  of 
Shockley-Read-Hall  recombination  and  Fowler-Nordheim 
tunneling.  Spontaneous  and  piezoelectric  polarizations  are 
modeled  as  a  quadratic  function  of  Al  mole  fraction  [10],  with 
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a  slight  linear  temperature  variation  or  pyroelectric  coefficient 
[10].  Acceptor- type  traps  at  a  uniform  volume  density  have 
been  placed  in  the  GaN  buffer. 

The  small-scale  simulation  includes  the  coupled  solution  of 
the  hydrodynamic  equations  for  electronic  transport  (a  general¬ 
ization  of  the  drift- diffusion  model)  and  the  thermal  diffusivity 
equation.  In  this  model,  heat  is  generated  as  electrons  thermal- 
ize  mostly  where  the  device  electric  field  is  typically  highest, 
i.e.,  in  the  channel  region  near  the  drain  edge  of  the  gate.  It 
was  found  that  this  made  a  clear  difference  for  simulations 
with  deep  traps  versus  drift-diffusion;  hot  electrons  diffuse  out 
of  the  channel  and  thermalize  into  the  deep  traps  in  the  GaN 
epilayer,  especially  for  the  highest  drain  biases  we  tested  here 
(65  V),  where  this  can  be  through  the  entire  layer  under  the 
channel  between  the  gate  and  drain.  More  sophisticated  models 
exist  such  as  Monte  Carlo  [11],  where  electrons  are  treated  as 
discrete  particles,  instead  of  an  electronic  fluid,  and  even  hot 
phonon  models  [12],  where  electron  energy  loss  is  modulated 
by  the  density,  energy  spectrum,  and  type  of  phonons  in  the 
crystal.  However,  this  requires  much  greater  computational 
complexity  and  thus  limits  the  physical  complexity  reasonably 
attainable.  Part  of  our  goal  was  to  build  a  versatile  model  based 
on  realistic  device  dimensions  at  all  relevant  length  scales. 

The  large-scale  model  solves  only  the  thermal  diffusivity 
equation  with  temperature-dependent  thermal  conductivities  for 
the  relevant  materials,  but  because  this  requires  a  3-D  solu¬ 
tion  of  a  multifinger  part,  there  is  considerable  computational 
and  complexity  savings  not  attempting  the  full  electrothermal 
solution. 

This  electrothermal  model  was  solved  in  2-D  for  portions  of 
the  device  including  the  source,  gate,  and  drain  metallizations; 
SiN  passivation  and  a  source-connected  field  plate,  the  epilay- 
ers;  and  about  10  fim  of  the  SiC  substrate,  in  order  to  include  all 
electrically  active  regions.  Thermal  boundary  conditions  at  the 
top  and  sides  are  adiabatic,  but  the  bottom  is  linked  to  a  thermal 
simulation  as  previously  described.  Specifically,  the  simulation 
domain  for  the  electrothermal  model  is  a  plane  perpendicular 
to  the  gate  finger,  which  is  20  /im  wide  at  the  epilayers,  with 
the  gate  at  the  center  and  including  some  of  the  source  and 
drain  ohmics.  This  domain  includes  a  portion  of  the  SiC  in 
the  unusual  shape  of  the  bottom  half  of  a  12-sided  regular 
polygon  with  the  SiC/epilayer  interface  directly  under  the  gate 
at  the  center.  The  reason  for  this  unusual  domain  is  driven 
by  the  physics  of  the  problem;  some  heat  is  being  generated 
in  the  channel,  but  most  is  at  the  drain  edge  of  the  gate.  At 
about  10  fi m  away  from  this  region  (which  is  a  distance  that  is 
considerably  greater  than  the  size  of  the  main  heat  source),  it 
is  observed  in  the  model  that  the  heat  is  radiating  away  such 
that  the  isothermal  contours  are  nearly  circular  and  centered 
about  the  drain  edge  of  the  gate.  With  at  most  3  °C  variation, 
the  temperature  of  this  unusual  SiC  bottom  domain  (chosen  to 
approximate  a  half-circle  but  with  fewer  mesh  elements)  was 
observed  to  be  constant  for  this  reason.  This  allows  a  very 
simple  connection  to  the  geometrically  complex  3-D  thermal 
model  to  be  described.  Specifically,  the  3-D  thermal  model 
was  run  for  different  power  levels  and  base-plate  temperatures 
(Tbp),  and  the  average  temperature  at  this  same  half-circle 
region  at  the  center  of  a  central  finger  of  the  multifinger  device 


Fig.  1.  (Thick  solid  lines)  Model  results  compared  with  (thin  dotted  lines) 
experimental  results  for  a  selected  device  from  wafer  “A”  with  better- than- 
average  pinch-off  characteristics.  Data  are  for  drain  bias  of  65,  48,  28,  and 
10  V  going  from  left  (dark  lines)  to  right  (light  lines). 


Fig.  2.  (a)  Impact  of  varying  the  piezoelectric  and  polarization  charge  at  the 

channel  on  device  performance,  (b)  Impact  of  varying  the  trap  density  in  the 
GaN  buffer.  For  both,  data  are  for  drain  bias  of  65, 48,  28,  and  10  V  going  from 
left  (dark  lines)  to  right  (light  lines). 


Vg(V) 


Fig.  3.  (Solid  curves)  Transconductance  curve  for  the  wafer  “A”  population 
of  devices  with  (thick  dashed  line)  overlaid  modeling  results. 


modeled  (Tnnk)  was  recorded.  The  center  was  chosen,  because 
most  of  the  device  is  expected  to  be  near  the  peak  temperature 
[2].  The  fit  of  the  temperature  at  this  half-circle  region  to 
a  thermal  resistance  as  a  function  of  only  Tavg  =  Tbp/2  + 
Tbnk/2  was  good  to  better  than  1.5  °C  error.  At  this  point, 
this  function  was  used  in  the  ISE  Sentaurus  Device  small-scale 
electrothermal  model  as  input  for  the  thermal  resistance  of  the 
bottom  domain.  The  average  temperature  along  this  region  of 
one  domain  agreed  with  the  average  temperature  along  the  same 
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TABLE  I 

Experiment  Versus  Model  Data  for  Several  Parameters 


Data  Set 

Id, max 
(raA/mm) 

1  D,SS 

(mA/mm) 

gmJVak 

(mS) 

gm.Peak 

Voltage  (V) 

Drain  Sensitivity 

A 

(V) 

Wafer  “A”  Data 
Median  ±  Sigma 

917  ±  56 

703  ±36 

261  ±  17 

-2.6  ±0.1 

0.1 17  ±  0.010 

-3.54  ±0.24 

Model  Data  for  “A” 

851 

672 

199 

-3.1 

0.128 

-3.91 

Process  “B”  Data 
Median  ±  Sigma 

930  ±41 

732  ±35 

224  ±  13 

-3.2  ±0.2 

0.129  ±  0.012 

-4.20  ±0.18 

Model  Data  for  “B” 

885 

731 

202 

-3.4 

0.144 

-4.32 

Wafer  “C”  Data 
Median  ±  Sigma 

642  ±  50 

323  ±31 

206  ±9 

-1.0  ±0.1 

0.105  ±  0.014 

-1.65  ±  0.14 

Model  Data  for  “C” 

657 

335 

208 

-0.92 

0.106 

-1.68 

Populations  of  devices  are  n=34  for  wafer  “A”  n=937  for  process  “B”,  and  n=98  for  wafer  “C’\  Because  the  data  is  non-normally 
distributed,  the  data  reported  is  the  median  and  the  errors  reported  are  calculated  using  the  interquartile  range  /  1 .35.  In  most  cases, 
this  gave  a  smaller  uncertainty  than  the  calculated  standard  deviation.  ID  max  is  at  gate  forward  bias  such  that  It;  =  ImA/mm  at 
Vo  =  10  V,  hyss  is  at  Vo  =  0  V  and  Vo  =  10  V,  and  drain  sensitivity  is  linear  slope  of  the  change  in  Vo  at  which  I(i  =  ImA/mm 
with  VD  over  the  voltage  range  [0,  15]. 


region  of  the  other  by  better  than  1  °C.  Because  the  models  are 
linked  in  this  way,  they  do  not  need  to  share  mesh  elements  as 
is  commonly  needed  in  multiscale  modeling. 

III.  Model  Calibration 

The  model  was  calibrated  using  two  wafers  that  were  run 
through  a  nearly  identical  process  with  one  process  variation 
between  them.  This  variation  was  expected  to  introduce  some 
difference  in  the  2-DEG  charge  density  under  the  gate,  due  to 
the  relationship  between  the  surface  potential  and  the  2-DEG 
charge  density,  but  otherwise  was  not  expected  to  impact  the 
epitaxy  stack,  structure,  and  cross  section  of  the  final  device 
geometry.  These  two  wafers  are  denoted  as  “A”  and  “B”  in  this 
paper.  As  a  measure  of  the  predictive  ability  of  this  model  and 
its  versatility  with  respect  to  its  applicability  to  different  epitax¬ 
ial  starting  materials,  the  simulated  electrical  characteristics  of 
a  third  wafer,  which  was  denoted  as  wafer  “C,”  are  also  shown. 
The  degree  of  agreement  between  the  simulated  and  measured 
electrical  characteristics  of  wafer  “C”  is  reported.  Wafer  “C” 
is  processed  in  identical  fashion  to  wafer  B;  however,  it  is 
different  in  its  epitaxial  stack,  which  is  described  in  Section  IV 
of  this  paper. 

Most  input  parameters  are  either  fairly  well  known  or  have 
only  a  modest  effect  on  the  electrical  properties  within  their 
expected  range  of  uncertainty.  However,  a  small  percentage 
change  in  the  polarization  charge  at  the  AlGaN/GaN  interface 
is  observed  to  have  a  large  effect  on  the  threshold  voltage  and 
saturation  current  of  the  modeled  device.  In  addition,  changes  in 
the  density  of  active  traps  in  the  GaN  buffer  affect  the  threshold 
voltage  and  saturation  current  as  before  but  also  affect  the 
subthreshold  slope  of  the  device.  This  is  of  some  importance 
from  a  model  perspective,  because  an  accurate  density  of  active 
traps  is  hard  to  quantitatively  determine  by  experiment.  In  both 
of  these  cases,  the  input  parameters  were  adjusted  to  allow 
more  exact  fitting  of  the  experimental  data.  For  the  polarization 
charge  term,  the  interface  charge  density  extracted  from  [10]  is 
multiplied  by  a  “charge  multiplier”  (CM)  close  to  1.0  to  fit  data, 
and  for  the  trap  density,  this  value  is  adjusted,  always  assuming 


that  the  trap  density  is  spatially  uniform  and  with  unchanged 
energy  level  and  capture  cross  section. 

Fig.  1  shows  subthreshold  data  collected  for  a  good  device 
on  wafer  “A”  at  several  voltages,  with  modeling  data  overlaid. 
Both  parameters  were  specifically  adjusted  to  match  the  se¬ 
lected  “better  performing”  device  as  well  as  possible,  meaning 
one  that  pinched  off  with  high  subthreshold  slope.  It  is  clear 
that  there  are  some  leakage  paths  that  are  not  included  in  the 
model.  A  better  fit  could  have  been  obtained  by  varying  the 
trap  density  with  depth  (traps  deep  in  the  buffer  and  far  from 
the  channel  were  seen  to  mainly  affect  only  the  higher  drain 
voltage  data),  but  this  was  not  done  so  as  to  minimize  adjustable 
parameters.  Fig.  2(a)  and  (b)  shows  the  effect  of  increasing 
the  CM  and  increasing  deep  levels,  respectively.  It  is  seen 
that  modest  changes  in  either  of  these  causes  big  changes  in 
the  electrical  properties  of  the  model  and  that  they  affect  the 
electrical  properties  in  very  different  way,  so  that  their  effects 
are  separable.  Fig.  3  shows  the  transconductance  curve  for  the 
same  wafer  “A”  devices  with  overlaid  modeling  results  for  the 
same  model  as  Fig.  1,  and  Table  I  shows  the  level  of  parametric 
data  agreement  seen,  which  were  generated  after  the  model 
calibration  was  complete.  Model  data  agree  within  two  standard 
deviations  for  all  data  but  those  related  to  gm,  Peak  for  wafer 
“A”.  As  can  be  seen  in  Fig.  3,  this  is  due  to  the  unusual  peak 
in  gm  for  the  experimental  data  for  this  wafer.  We  are  unaware 
of  the  precise  cause  of  this  anomaly,  which  was  not  present  on 
some  other  wafers. 


IV.  Extended  Model  Calibration 

Fig.  4  shows  data  from  wafer  “B”  run  with  the  same  device 
geometry  but  under  a  different  process;  it  is  experimentally  seen 
that  the  data  appear  to  form  a  different  data  set  and  to  require 
a  new  model.  The  solid  line  represents  data  from  the  model 
after  changing  only  the  CM  and  illustrates  the  importance 
of  controlling  process  parameters  that  affect  the  surface  of 
GaN  devices.  It  has  been  observed  before  that  the  electrical 
properties  of  GaN  devices  are  sensitive  to  the  process,  bias 
history,  and  past  light  exposure,  with  persistent  effects  seen  that 
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Vg(V) 


Fig.  4.  (Solid  curves)  Transconductance  curve  for  the  wafer  “B”  population 
of  devices  with  (thick  dashed  line)  overlaid  modeling  results. 


Vg(V| 


Fig.  5.  (Solid  curves)  Transconductance  curve  for  the  wafer  “C”  population  of 
devices  with  thinner  barrier  and  lower  mole  fraction  AlGaN  with  (thick  dashed 
line)  overlaid  modeling  results. 

are  most  commonly  attributed  to  energetically  deep  traps  that 
affect  the  density  of  charge  in  the  channel  [13],  [14].  As  could 
be  expected,  because  the  transconductance  data  of  Fig.  4  are 
less  noisy  than  those  of  Fig.  3,  the  parametric  data  in  Table  I 
agree  better  with  the  model. 

At  this  point,  the  proposed  model  can  be  used  for  predictive 
purposes.  The  AlGaN  barrier  in  the  model,  as  matched  to 
wafer  “B,”  was  changed  to  correspond  to  an  experimental  wafer 
“C”  with  reduced  AlGaN  barrier  thickness  and  mole  fraction. 
Nothing  else  was  changed,  including  the  CM  in  this  instance. 
Fig.  5  shows  the  agreement  between  model  and  experiment,  and 
Table  I  shows  the  parametric  agreement. 

V.  Conclusion 

A  physics-based  finite-element  model  containing  the  entire 
electrical  region  of  interest  for  dc  operation  and  the  entire 
thermal  region  including  the  package  has  been  generated.  The 
model  has  been  calibrated  to  several  dc  metrics,  and  the  level  of 
agreement  has  been  shown.  It  has  been  shown  that  this  model 
can  be  easily  adjusted  to  fit  additional  populations  of  devices, 
keeping  the  fundamental  physics  in  the  model  the  same.  It  has 
also  been  shown  that  the  electrical  properties  of  the  model  are 
sensitive  to  trap  density  within  the  GaN  buffer  layer,  helping 
constrain  this  parameter. 
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